Hyperparameters and Model Validation

Previously, we saw the basic recipe for applying a supervised machine learning model:

Choose a class of model
Choose model hyperparameters
Fit the model to the training data
Use the model to predict labels for new data

The first two pieces of this are the most important part of using these tools and techniques effectively.
The question that comes up is: How can we make informed choice for these parameters?
We've touched upon questions from this realm already, but here we are going to examine it in more detail.

Thinking about Model Validation

In principle, model validation is very simple:

choosing a model and its hyperparameters
Estimation: applying it to some of the training data and comparing the prediction to the known value.

Let's first go through a naive approach and see why it fails.

Model validation the wrong way

Let's demonstrate the naive approach to validation using the Iris data, which we saw previously:



In [ ]:

    
from sklearn.datasets import load_iris
iris = load_iris()
X = iris.data
y = iris.target

Next we choose a model and hyperparameters:



In [ ]:

    
from sklearn.neighbors import KNeighborsClassifier
model = KNeighborsClassifier(n_neighbors=1)

Then we train the model, and use it to predict labels for data we already know:



In [ ]:

    
model.fit(X, y)
y_model = model.predict(X)

Finally, we compute the fraction of correctly labeled points:



In [ ]:

    
from sklearn.metrics import accuracy_score
accuracy_score(y, y_model)

We see an accuracy score of 1.0, which indicates that 100% of points were correctly labeled by our model!
But is this truly measuring the expected accuracy?
Have we really come upon a model that we expect to be correct 100% of the time?

Model validation the right way: Holdout sets

We need to use a holdout set: that is, we hold back some subset of the data from the training of the model, and then use this holdout set to check the model performance.
This splitting can be done using the train_test_split utility in Scikit-Learn:



In [ ]:

    
from sklearn.model_selection import train_test_split
# split the data with 50% in each set
X1, X2, y1, y2 = train_test_split(X, y, random_state=0, train_size=0.5)

# fit the model on one set of data
model.fit(X1, y1)

# evaluate the model on the second set of data
y2_model = model.predict(X2)
accuracy_score(y2, y2_model)

The nearest-neighbor classifier is about 90% accurate on this hold-out set, which is more inline with out expactation.
The hold-out set is similar to unknown data, because the model has not "seen" it before.

But, we have lost a portion of our data to the model training.
In the above case, half the dataset does not contribute to the training of the model! This is not optimal.

Model validation via cross-validation

One way to address this is to use cross-validation; that is, to do a sequence of fits where each subset of the data is used both as a training set and as a validation set. Visually, it might look something like this:

Here we do two validation trials, alternately using each half of the data as a holdout set. Using the split data from before, we could implement it like this:



In [ ]:

    
y2_model = model.fit(X1, y1).predict(X2)
y1_model = model.fit(X2, y2).predict(X1)
accuracy_score(y1, y1_model), accuracy_score(y2, y2_model)

We could compute the mean of the two accuracy scores to get a better measure of the global model performance.
This particular form of cross-validation is a two-fold cross-validation.

We could expand on this idea to use even more trials, and more folds.
Here is a visual depiction of five-fold cross-validation:

We split the data into five groups
Each of them is used to evaluate the model fit on the other 4/5 of the data.
What is the advantage of higher-degree crossvalidation? What is the drawback?

This would be rather tedious to do by hand, and so we can use Scikit-Learn's cross_val_score convenience routine to do it succinctly:



In [ ]:

    
from sklearn.model_selection import cross_val_score
cross_val_score(model, X, y, cv=5)

This gives us an even better idea of the performance of the algorithm.

How could we take this idea to its extreme?

The case in which our number of folds is equal to the number of data points.
This type of cross-validation is known as leave-one-out cross validation, and can be used as follows:



In [ ]:

    
from sklearn.model_selection import LeaveOneOut
import numpy as np

loo = LeaveOneOut()
loo.get_n_splits(X)
scores = []
for train_index, test_index in loo.split(X):
    model.fit(X[train_index], y[train_index])
    scores.append(accuracy_score(y[test_index], model.predict(X[test_index])))
scores = np.array(scores)
scores

Because we have 150 samples, the leave one out cross-validation yields scores for 150 trials, and the score indicates either successful (1.0) or unsuccessful (0.0) prediction. Taking the mean of these gives an estimate of the error rate:



In [ ]:

    
scores.mean()

This gives us a good impression on the performance of our model. But there is also a problem. Can you spot it?

Selecting the Best Model

Now that we can evaluate a model's performance, we will tackle the question how to select a model and its hyperparameters.
A very important question to ask is: If my estimator is underperforming, how should I move forward? What are the options?

Use a more complicated/more flexible model
Use a less complicated/less flexible model
Gather more training samples
Gather more data to add features to each sample

Sometimes the results are counter-intuitive:

Using a more complicated model will give worse results
Adding more training samples may not improve your results

The Bias-variance trade-off

Fundamentally, the question of "the best model" is about finding a sweet spot in the tradeoff between bias and variance. Consider the following figure, which presents two regression fits to the same dataset:

Comment on the models.
Which model is better?
Is either of the models 'good'? Why?

Consider what happens if we use these two models to predict the y-value for some new data.
In the following diagrams, the red/lighter points indicate data that is omitted from the training set:

It is clear that neither of these models is a particularly good fit to the data, but they fail in different ways.

The score here is the $R^2$ score, or coefficient of determination, which measures how well a model performs relative to a simple mean of the target values.

$R^2=1$ indicates a perfect match
$R^2=0$ indicates the model does no better than simply taking the mean of the data
Negative values mean even worse models.

Left model

Attempts to find a straight-line fit through the data.
The data are intrinsically more complicated than a straight line, the straight-line model will never be able to describe this dataset well.
Such a model is said to underfit the data: that is, it does not have enough model flexibility to suitably account for all the features in the data
Another way of saying this is that the model has high bias

Right model

Attempts to fit a high-order polynomial through the data.
The model fit has enough flexibility to nearly perfectly account for the fine features in the data
It very accurately describes the training data
Its precise form seems to be more reflective of the noise properties rather than the intrinsic properties of whatever process generated that data
Such a model is said to overfit the data: that is, it has so much model flexibility that the model ends up accounting for random errors as well as the underlying data distribution
Another way of saying this is that the model has high variance.

From the scores associated with these two models, we can make an observation that holds more generally:

For high-bias models, the performance of the model on the validation set is similar to the performance on the training set, but the overall score is low.
For high-variance models, the performance of the model on the validation set is far worse than the performance on the training set.

Imagine we have the ability to tune the model complexity, then we can expect the training score and validation score to behave as illustrated in the following figure:

The diagram shown here is often called a validation curve, and we see the following essential features:

The training score is everywhere higher than the validation score. This is generally the case. Why?
For very low model complexity (a high-bias model), the training data is under-fit, which means that the model is a poor predictor both for the training data and for any previously unseen data.
For very high model complexity (a high-variance model), the training data is over-fit, which means that the model predicts the training data very well, but fails for any previously unseen data.
Which model/level of complexity should we choose?.

Validation curves in Scikit-Learn

Let's look at an example. We will use a polynomial regression model: this is a generalized linear model in which the degree of the polynomial is a tunable parameter. For example, a degree-1 polynomial fits a straight line to the data; for model parameters $a$ and $b$:

$$ y = ax + b $$

A degree-3 polynomial fits a cubic curve to the data; for model parameters $a, b, c, d$:

$$ y = ax^3 + bx^2 + cx + d $$

We can generalize this to any number of polynomial features.
In Scikit-Learn, we can implement this with a simple linear regression combined with the polynomial preprocessor.
We will use a pipeline to string these operations together:



In [ ]:

    
from sklearn.preprocessing import PolynomialFeatures
from sklearn.linear_model import LinearRegression
from sklearn.pipeline import make_pipeline

def PolynomialRegression(degree=2, **kwargs):
    return make_pipeline(PolynomialFeatures(degree), LinearRegression(**kwargs))

Now let's create some data to which we will fit our model:



In [ ]:

    
import numpy as np

def make_data(N, err=1.0, rseed=1):
    # randomly sample the data
    rng = np.random.RandomState(rseed)
    X = rng.rand(N, 1) ** 2
    y = 10 - 1. / (X.ravel() + 0.1)
    if err > 0:
        y += err * rng.randn(N)
    return X, y

X, y = make_data(40)

We can now visualize our data, along with polynomial fits of several degrees:



In [ ]:

    
%matplotlib inline
import matplotlib.pyplot as plt
import seaborn; seaborn.set()  # for beautiful plotting

X_test = np.linspace(-0.1, 1.1, 500)[:, None]

plt.scatter(X.ravel(), y, color='black') # plot data
axis = plt.axis()

# plot ploynomials
for degree in [1, 3, 5]:
    y_test = PolynomialRegression(degree).fit(X, y).predict(X_test)
    plt.plot(X_test.ravel(), y_test, label='degree={0}'.format(degree))

plt.xlim(-0.1, 1.0)
plt.ylim(-2, 12)
plt.legend(loc='best');

We can controll the model complexity (the degree of the polynomial), which can be any non-negative integer.
A useful question to answer is this: what degree of polynomial provides a suitable trade-off between bias (under-fitting) and variance (over-fitting)?
We can make progress in this by visualizing the validation curve
This can be done straightforwardly using the validation_curve convenience routine
Given a model, data, parameter name, and a range to explore, this function will automatically compute both the training score and validation score across the range:



In [ ]:

    
from sklearn.model_selection import validation_curve
degree = np.arange(0, 21)
train_score, val_score = validation_curve(PolynomialRegression(), X, y,
                                          'polynomialfeatures__degree', degree, cv=10)

plt.plot(degree, np.median(train_score, 1), color='blue', label='training score')
plt.plot(degree, np.median(val_score, 1), color='red', label='validation score')
plt.legend(loc='best')
plt.ylim(0, 1)
plt.xlabel('degree')
plt.ylabel('score');

This shows precisely the behavior we expect:

The training score is everywhere higher than the validation score
The training score is monotonically improving with increased model complexity
The validation score reaches a maximum before dropping off as the model becomes over-fit.

Which degree ist optimal?

Let's plot this.



In [ ]:

    
plt.scatter(X.ravel(), y)
lim = plt.axis()
y_test = PolynomialRegression(3).fit(X, y).predict(X_test)
plt.plot(X_test.ravel(), y_test);
plt.axis(lim);

Learning Curves

One important aspect of model complexity is that the optimal model will generally depend on the size of your training data.
For example, let's generate a new five times larger dataset:



In [ ]:

    
X2, y2 = make_data(200)
plt.scatter(X2.ravel(), y2);

We will duplicate the preceding code to plot the validation curve for this larger dataset; for reference let's over-plot the previous results as well:



In [ ]:

    
degree = np.arange(51)
train_score2, val_score2 = validation_curve(PolynomialRegression(), X2, y2,
                                            'polynomialfeatures__degree', degree, cv=10)

plt.plot(degree, np.median(train_score2, 1), color='blue', label='training score')
plt.plot(degree, np.median(val_score2, 1), color='red', label='validation score')
plt.plot(degree[:train_score.shape[0]], np.median(train_score, 1), color='blue', alpha=0.3, linestyle='dashed')
plt.plot(degree[:train_score.shape[0]], np.median(val_score, 1), color='red', alpha=0.3, linestyle='dashed')
plt.legend(loc='best')
plt.ylim(0, 1)
plt.xlabel('degree')
plt.ylabel('score');

From the validation curve it is clear that the larger dataset can support a much more complicated model:

The peak here is probably around a degree of 6
Even a degree-20 model is not seriously overfitting the data

Thus we see that the behavior of the validation curve has not one but two important inputs:

the model complexity.
the number of training points.

It is often useful to explore the behavior of the model as a function of the number of training points.
This can be done by using increasingly larger subsets of the data to fit our model.
This is called a learning curve.

The general behavior of learning curves is this:

A model of a given complexity will overfit a small dataset: this means the training score will be relatively high, while the validation score will be relatively low.
A model of a given complexity will underfit a large dataset: this means that the training score will decrease, but the validation score will increase.
A model will never, except by chance, give a better score to the validation set than the training set: this means the curves should keep getting closer together but never cross.

With these features in mind, we would expect a learning curve to look qualitatively like that shown in the following figure:

The Learning curve converges to a particular score as the number of training samples grows
Once we have enough training data a particular model has converged
That means: adding more training data will not help!
To increase the model performance another (often more complex) model must be chosen

Learning curves in Scikit-Learn

Scikit-Learn offers a convenient utility for computing such learning curves:



In [ ]:

    
from sklearn.model_selection import learning_curve

fig, ax = plt.subplots(1, 2, figsize=(16, 6))
fig.subplots_adjust(left=0.0625, right=0.95, wspace=0.1)

for i, degree in enumerate([2, 9]):
    N, train_lc, val_lc = learning_curve(PolynomialRegression(degree),
                                         X, y, cv=10,
                                         train_sizes=np.linspace(0.3, 1, 25))

    ax[i].plot(N, np.mean(train_lc, 1), color='blue', label='training score')
    ax[i].plot(N, np.mean(val_lc, 1), color='red', label='validation score')
    ax[i].hlines(np.mean([train_lc[-1], val_lc[-1]]), N[0], N[-1],
                 color='gray', linestyle='dashed')

    ax[i].set_ylim(0, 1)
    ax[i].set_xlim(N[0], N[-1])
    ax[i].set_xlabel('training size')
    ax[i].set_ylabel('score')
    ax[i].set_title('degree = {0}'.format(degree), size=14)
    ax[i].legend(loc='best')

This plot gives us a visual depiction of how our model responds to increasing training data.
When the learning curve has already converged (i.e., when the training and validation curves are already close to each other) adding more training data will not significantly improve the fit!.
This situation is seen in the left panel, with the learning curve for the degree-2 model.
To increase the converged score a more complicated model must be used.
In the right panel: by moving to a much more complicated model, we increase the score of convergence
The drawback is a higher model variance (indicated by the difference between the training and validation scores).
If we were to add even more data points, the learning curve for the more complicated model would eventually converge.

Plotting a learning curve for your particular choice of model and dataset can help you to make this type of decision about how to move forward in improving your analysis.

What is the difference between a validation curve and a learning curve?

Validation in Practice: Grid Search

The preceding discussion gave you some intuition into the trade-off between bias and variance, and its dependence on model complexity and training set size.
In practice, models generally have more than one knob to turn, and thus plots of validation and learning curves change from lines to multi-dimensional surfaces.
Such visualizations are difficult.
We would like to find the particular model that maximizes the validation score directly, instead.

Scikit-Learn provides automated tools to do this in the grid search module.
Here is an example of using grid search to find the optimal polynomial model.
We will explore a three-dimensional grid of model features:

Polynomial degree
The flag telling us whether to fit the intercept
The flag telling us whether to normalize the problem.

This can be set up using Scikit-Learn's GridSearchCV meta-estimator:



In [ ]:

    
from sklearn.model_selection import GridSearchCV

param_grid = {'polynomialfeatures__degree': np.arange(21),
              'linearregression__fit_intercept': [True, False],
              'linearregression__normalize': [True, False]}

grid = GridSearchCV(PolynomialRegression(), param_grid, cv=10)

Notice that like a normal estimator, this has not yet been applied to any data.
Calling the fit() method will fit the model at each grid point, keeping track of the scores along the way:



In [ ]:

    
grid.fit(X, y);

Now that this is fit, we can ask for the best parameters as follows:



In [ ]:

    
grid.best_params_

Finally, if we wish, we can use the best model and show the fit to our data using code from before:



In [ ]:

    
model = grid.best_estimator_

plt.scatter(X.ravel(), y)
lim = plt.axis()
y_test = model.fit(X, y).predict(X_test)
plt.plot(X_test.ravel(), y_test);
plt.axis(lim);